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Abstract 

We name strongly localized quantum crystalline state (SLQCS) a de- 
terminantal wave-function of the single-particle wave-functions obtained 
by crystalline translations of a wave-function different from zero only 
within a primitive cell of the considered quantum crystalline phase. SLQCSs 
accurately reproduce the low density behaviour of the quantum crystals 
that, as Wigner's crystalline phase of the jellium model, can become 
very dilute. Our analysis explicitly deals with this system. We show 
that the SLQCS energy per particle at large dilution (r^) behaves as 
—Mdi,a /ts + Ca/rs^^'^ , where Mdi,cr turns out to be the Madelung con- 
stant of the considered cubic symmetry a and Ca is a positive constant 
only numerically determined. Moreover, as the density gets smaller and 
smaller, each electron becomes more and more confined to the centre of 
its cell. 
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1 Introduction 



The jellium model [1] [2] is the simplest model of metallic conductors. It ap- 
proximates the crystalline array of positive ions by a classical, uniform, positive 
charge density while the valence electrons are described as quantum mechan- 
ical particles interacting among themselves and with the electric field of the 
aforesaid neutralizing charge distribution. In this way, the Hamiltonian of the 
system, confined at first within a box of volume V , is 

2 2 I- I- -^l\r-T'\ JX, r p-/i|ri-r| 

2 Jy Jy |r-r'| ^Jy \r, - v] 

_V^ + - V 

^ 2m 2 ^ |ri -r,| ' 

i=l l<i^j<N ' 

Here, N denotes the number of the electrons contained in the box, n = N/V 
the electron number density, the position vector of the ith electron, e and m 
electrons' charge and mass. The cut-off parameter fj, goes to zero after taking 
the thermodynamic limit: V — > cx), N ^ oo with n fixed. Besides, in Eq. ([1]) 
and throughout the paper we adopt units such that h = 1. 

The first order perturbative approximation of the energy of the fundamental 
state of H was evaluated by Wigner^Sj for the unpolarized case and reads 

where denotes the fundamental state of N free and unpolarized spin 1/2 
particles, ag = 1/me^ is the Bohr atomic radius, = ro/ao the perturbation 
parameter with rg = (3/47rn)^/'^. The corresponding approximation for the 
polarized case was obtained by Bloch[4 and reads 

where \Fp) denotes the fundamental state of N polarized free fermions. Expres- 
sions ^ and ^ refer to homogeneous systems. Their comparison shows that 
the polarized configuration is the one stable at lower density, i.e. for > 5.7. 
The energy of the classical crystals, formed by electrons arranged on a lattice 
with simple (sc) or body (bcc) or face centred cubic (fee) symmetry and a 
neutralizing positive uniform background, is —N(e'^ /2ao)Mdi^a-/'rs, where Mdi^a 
with <7 — 1,2,3 denotes the Madelung constant respectively equal to 1.76012, 
1.79186 or 1.79175 5 for the sc, bcc or fee case. Since at high dilution these 
energies become definitely smaller than Eq. ([3])'s, Wigner concluded that, at 
large Ts, the fundamental state of the jellium no longer can be that of the ho- 
mogeneous Fermi gas at T = {f'K. It must have a crystalline structure. In 
fact, assuming that each electron oscillates around its equilibrium position and 



2 



that these positions form a bee lattice, Wigner[B] obtained the following energy 
expression 



where the term proportional to rj^, representing the energy of the classical 
crystal, was put in, so to say, by hand. A better approximation was later 
obtained by Carr[7J expanding, as Wigner did, the electrons' position vectors 
around the bcc lattice values and, contrarily to Wigner, also retaining the lowest 
order terms that couple the resulting harmonic oscillators among themselves. 
The resulting first order approximation of the fundamental state energy is 



cannot however invoke the Ritz-Rayleigh principle to conclude that at high di- 
lution the jellium model is in a bcc crystalline state because expressions Q and 
([5|) were not obtained considering the full Hamiltonian H of the jellium model, 
as defined in Eq. (1.1), but two different modifications of it. The property can 
however be proved showing that the solution, with the bcc symmetry, of the 
jellium's Hartree-Fock (HF) equation yields an energy lower than Eq. (I3])'s at 
large r-g. To this aim we recall that the solution ■(/'(ri, . . . , rjv) of the HF equa- 
tions is the wave function that has a determinantal structure and minimizes 
the expectation value of if[51[S]. Further, in writing down the determinantal 
wave-function, one usually assumes that the single particle wave-functions have 
a particular functional form. In particular, assuming that the single particle 
wave- functions are the plane waves e^^^ '^ with j = 1, . . . , A^, the resulting HF 
energy coincides with expressions ([2]) or ([3]) depending on the considered po- 
larization degree. On the contrary, a HF crystalline solution can be obtained 
assuming that the single particle wave-functions ^Pjiv) have the form e*'^J ''u(r), 
with j = 1, . . . , TV and u{r) a periodic function (with the chosen crystalline sym- 
metry) to be determined through the minimization of the HF energy. Applying 
this procedure, it was recently found^lOj that crystalline HF energies become 
smaller than Eq. Q's at > 4.4 and that the bcc HF energy is the smallest 
in the region > 13.3. In this way, Wigner's statement that diluted jellium 
assumes a bcc crystalline configuration is put on a firmer basis because it is 
now a consequence of the Ritz-Rayleigh principle. Strictly speaking the stabil- 
ity of the crystalline phase cannot be considered rigorously proved because the 
exact fundamental state is not required to have a determinantal structure, as 
assumed in the HF equations. However, the statement's validity is also sup- 
ported by Quantum Monte Carlo (QMC) calculations Hll fT2 l [151 [Ml [T5 l flB] 
that, relaxing the determinantal restriction, confirm that the stable phase of 
the dilute electron jellium is the bcc one (even though the exact value of the 
transition is not yet well assessed because, depending on the paper, it ranges 
from 50 to 105). 







At large we have EcrC < E^rW < Ep. On the basis of these inequalities one 
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The aim of this paper is to present a simpler procedure to show the existence, 
for the jelhum model, of crystalline phases more stable than the fluid one and 
to illustrate how these quantum-mechanical crystalline solutions reproduce the 
classical crystal behaviours as — > oo. These solutions are essentially obtained 
by a simplified HF procedure as follows. First, one assumes that the confining 
box V consists of = {2M + 1)^ primitive cells with the chosen crystalline 
symmetry. Then, one assumes that the single particle wave-function ipj{r,a), 
with j = 1, . . . , N, differs from zero only within the jth primitive cell and goes 
to zero at the cell's border. Besides, each (pj{r,a) is obtained translating a 
normalized function (j)[r, a) that differs from zero only within Vq (the primitive 
cell centred at the origin) and goes to zero at the border of Vq. The correspond- 
ing translation vector is the a^s being the vectors specifying the 
primitive cell and nij = f^j.s) the vector (with integer components) 
that labels the jth cell. Finally, a is a parameter to be determined by a varia- 
tional procedure that will be described later. It must be noted that, if j ^ k, 
functions a) and </5fc(r, a) are orthogonal because their supports have void 
intersection. The determinantal wave-function is 

■(/'(ri,...,rAr,a) = — L= ^(-l)^(^i(r,, , a)(y92(ri2 , a) . . . (^Ar(ri„ , a), (6) 
V jV. p 

where the sum is performed over all the permutations ii, Z2, . . . , zat of 1, 2, . . . , TV 
and (—1)^ is the parity of the considered permutation. Hereafter, wave-function 
© is referred to as a strongly localized quantum crystalline state (SLQCS) if the 
involved ipj(r^a)s obey the properties mentioned above. Then, one evaluates 
the expectation value of H over the above normalized wave- function. In the 
thermodynamic and subsequent /i — > limits, the resulting energy expression 
depends on a and r^, and one looks for the a value that makes the energy 
value minimum for each value. In this way one finds (see Fig.l) that: i) 
the resulting sc, bcc and fee SLQCS energies become smaller than Eq. ^'s at 
Ts > 28, Vs > 36 and > 38, respectively, ii) the SLQCS fee energy nearly 
coincides with the bcc one. More definitely, it is always smaller than the bee's 
and their slight difference monotonically decreases throughout the explored 
range [0, 10*]. The sc SLQCS energy at first is smaller that the fee's and bee's 
and becomes greater of the last two at > 500, iii) the three SLQCS energies 
approach the relevant Madelung contributions as — > oo. Consequently, the 
bcc SLQC energy approaches the corresponding leading term of Eq.s (|4|) and 

iv) the optimized a value increases with and at large one has a oc 

v) in the same region, the subleading asymptotic terms of the sc, bcc and 
fee SLQCS energies are proportional to rj^^^, and vi) our SLQCS energies are 
always greater than those obtained in Ref.[TD] that lie close to the full circles 
shown in Fig. la. 

Result vi) is not unexpected because the above HF procedure minimizes 
{tjj\H\4') over a class of functions that, for the assumption that the (pj{r,a)s 
have compact support, is more restricted than that used in Ref. 10 . The SLQCS 
procedure is however interesting for its numerical simplicity and for the choice of 
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Figure 1: Left and right parts: respectively the energy per particle of the fully 
polarized jellium model (in Rydberg units) in the inner and outer part of the 
Ts-range. (Note the logarithmic-scale on the right.) In both figures, the 
continuous hues refer to Eq. the dotted hues to Carr's bcc expression ([5]), 
and the short-dash and the long-dash curves to the SLQC solutions with the sc 
or the bcc symmetry, respectively. The full triangles and circles correspond to 
the QMC results of Ref.s [HI [16]. The SLQCS fee solution is not shown because 
it lies very close to the bcc. 
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the single particle wave-functions. The latter's noticed orthogonality property 
has an important consequence: if we evaluate the expectation value of H over 
the normalized Hartree wave-function 

V'/f(ri, ...,TN,a)^ (pi(ri,a)v32(r2,a) • • ■ </?Ar(rAr, a), (7) 

we obtain the same expectation value resulting from the fully antisymmetric 
state (O. In other words, SLQCSs have the property that the Hartree and the 
Hartree-Fock equations coincide because the exchange contributions turn out to 
be equal to zero as it will be explicitly shown below Eq. ([^5]) . The equivalence of 
Eq. ([6]) with ([7]) in their final results and the condition that the '/3j(r, q;)s vanish 
on their cell borders imply that in the the SLQCS approach each electron in 
practice is bound to lie within a single cell of the lattice. This feature on the 
one hand explains why wave-function ^ was named strongly localizedquantum 
crystalline state. On the other hand, the same feature is shown by classical ionic 
crystal. Hence, it is not surprising that the quantum-mechanical description, 
based on SLCQSs, reproduces the classical crystal behaviour in the limit — s- 
00. This appears already evident from results iii) and iv). In particular result 
iv) implies that, as rs — > oo, each electron becomes more and more confined 
to the centre of the primitive cell that the electron occupies, while result iii) 
quantum mechanically reproduces the 0(l/rs) term, manually introduced into 
Eq.s (1.4) and (1.5) by Wigner and Carr. 

The derivation of these results will now be detailed according to the following 
plan. In section 2 we work out the analytical expressions of the expectation 
value {il)\H\il>) in direct and reciprocal space for the three mentioned cubic 
symmetries. In section 3 we specify the functional form of 0(r, a) and report 
the relevant numerical results already anticipated at i)-vi). Section 4 concludes 
the paper. Some mathematical details are left to appendices A and B. 



2 Analytical expressions 

The primitive cell Vq centred at the origin is defined as 
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l-o^{r r = ^^, l<e,<l, i = 1,2,3}, (8) 
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where the involved a^s are respectively specified in Eq.s (27a), (29) and (34) of 
Chapt. 2 of Ref. [17] for the sc, bcc and fee lattices. For these symmetries the 
primitive cell volumes are 

vo^a\ aV2, aV4, (9) 

a denoting the cubic lattice spacing. The corresponding values are 

Ts = (3/47r)i/^a/ao, {Z/Snf/^a/ aa and (3/167r)i/^a/ao. (10) 

We define now the wave-function relevant to a SLQCS. To this aim, we first 
introduce a real even function (/)(r, a) that differs from zero only within Vb 
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and depends on a positive real parameter a. Whatever a, </)(r, a) and its first 
and second partial derivatives are continuous throughout Vb and vanish as r 
approaches the border of Vq. Besides, 0(r, a) is normalized, i.e. 

f ^^ir,a)dv = l, (11) 

and the a dependence is such that 

lim 0^(r,a) = ^(r), (12) 

a — >oo 

S{r) denoting the three-dimensional (3D) Dirac function. Assume now that the 
box of volume where the jellium is confined, contains (2M + 1) primitive 
cells along each direction a^. Then the jellium will contain N — {2M + 1)^ 
electrons, Vq is one of the possible cells and vq = V/N. We denote the primitive 
cells either as Vm (with —M < rrii < M and i = 1,2, 3) or, after performing a 
convenient relabeling m j(ni), as V, with j — 1, . . . ,N. The inverse relabeling 
exists and reads j ni(j). We set now 

3 

(Pinir,a) = (Pjir,a) = (f){r - ^■miai,a). (13) 

1=1 

Function (fimi^, Oi) results from the translation of (l){r) by It differs 

from zero only within primitive cell Vm and is there normalized. Moreover the 
(Pmlr, a)s are orthonormal, i.e. 

(Pj{r,a)(pj,{r,a)dv ^ Sj.j,, j,f = l,...,N. (14) 

Jv 

A strongly localized quantum crystalline (SLQC) wave-function is defined as the 
determinantal function defined by Eq.s and We assume that the 

normalized wave function of the N polarized electrons is a SLQC wave-function 
that, by the previous definitions, is normalized within V. With such a function, 
the electron number density n(r) turns out to be periodic and equal to 

3 

n(r) = ^0^(r-f-^mja,,a). (15) 

m 1=1 

As already anticipated, our task now is to evaluate the expectation value E{N, V, a) 
of H over a SLQC wave-function and to choose a so as to make E{N, V, a) as 
small as possible for each value. After putting 

eirs,a) = E{N,V,a)/N = {'4j\H\,P)/N, (16) 

from Eq.s ^ and (fT^ one gets 

eirs,a) = hb^b + K-b + {f) + {V) (17) 



7 



with 

= - y E / ,^^3 jr-—^ (19) 

(f ) - / V0(r, a) • V0(r, a)dz;, (20) 

and 

(^) = ^ E (21) 

where [H] 

J,,j = / dt. / dv'^^\v,a)^j''{v\a f (22) 
■/y Jy |r — r I 

Xj.j EE / / duV2(i"7a)'^i(i"':a)<^j(i", a)(^j(r',a)- — . (23) 

Jv Jv |r — r I 

Contributions and /ig-b respectively originate from the first and second 

term on the right hand side (rhs) of Eq. ((T|). Property (fO]) has been used to 
derive Eq.s and (^0)) . Quantities Jjj and i^jj represent the direct and the 
exchange contributions, respectively. The exchange contributions originate from 
the antisymmetry of the wave-function. In our case it turns out that the i^ijS 
with I 7^ J, as the sum present in Eq. (I2ip requires, are equal to zero owing to 
the fact that fi{r, a) and (pj(r, a) have supports with a void intersection. Thus, 
despite the fact that the considered wave-function is antisymmetric, SLQCS are 
characterized by vanishing exchange contributions. In other words, SLQCSs are 
defined by the property that the associated Hartree and Hartree-Fock equations 
yield the same final results. In fact, starting from Eq. (1.7), one easily checks 
that the resulting expectation of H again is the sum of analytical expressions 
(HOD and (EZl). 

Using Eq.s and ([H]), the sum involving the J^^s becomes 

With the change (m; — m'l) ^ mi and letting N go to infinity one gets 

g-/i|r-r' + X;, = i miSLi 



V Aj-^nV]' [ dv [ dwV^(r,a)0^(r',a) 

, ' ™ JVo JVo 



l<i^j<N 



(24) 
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where m runs over all the points of the lattice excluding the origin. The last 
restriction is signaled by the prime. Recalling that V 
he-b can be converted into a double sum, i.e. 



I 1^ V 



contribution 



e 

v 



E 



dv / dv 

Vo Jvo 



,(/>2(r,a)e~'^l'-'"'+2^ti(™.-™'0ad 



m'i)ai\ 



In the thermodynamic limit one finds that 



E 



dv I dv 

Vo JVo 



r — r 



In a similar way one shows that 



b-b 



2 ? 



2 ^ 



Vo 



Vo 



|r-r' + X;Li™ia«l 



(25) 



(26) 



Collecting the above results and taking the limit fi ^ oo one obtains the contri- 
bution of the total Coulombic interaction to the quantum energy per particle, 
namely 



hb-b + he-b + {V) = 



dv / dv 

V„ JVr 



, 27i0^(r, a) — -n? 
r- r'l 

2/^ „,\ I „2 



dv I dv 

Vo Jv 



, (r, q;)(/) (r', a) — 2n<j)^{r, a) + n 
|r - r' + 



(27) 



We make now three remarks. First, if one assumes that (j)'^{r,a) is a Gaussian 
function, the rhs of Eq. ([27]) coincides with the expression reported by Ewald 
in his calculation of the Coulombic energy per particle for a classical crystal 
(see e.g. Ref.s [El [20]). However, the introduction of the Gaussian function 
in Ewald's procedure was only a trick to make the convergence faster, while in 
the SLQCS procedure 0^(r, a) is the quantum probability density of finding an 
electron at position r. Second, the series present in Eq. ((27)) is convergent. In 
fact, setting A(m)z>(m) = J2i=i ^i^i where ^'(m) is a unit vector, the expansion 
of the denominator, present in Eq. ((27l) . at large A{m) values yields 



1 



1 



|r — r' + yl(m)i>(m) 



A{m) 



1 - 



(r ~ r') • i>(m) 
A(m) 



3(r - rQ ■ ;)(m) - (r - rQ^ 
2^2 (m) 



(28) 



The volume integrals of the first term on the rhs is equal to zero. The remaining 
two terms also yield vanishing contributions provided the sum over m is first 
performed over all the ms that have different directions and equal modulus and, 
subsequently, over the different |m| values. In this way, as ^(m) becomes very 
large, the first of the resulting sums amounts to performing an angular average 



9 



over the directions of j>(m). It is easily checked that the i> dependence of the 
numerators present in Eq. ((28|) is such that their angular averages vanish for 
the sc, fee and bcc cases. One concludes that the terms of series (l?7|) decrease 
as A~'^{m.) at large ms, and the series convergence is proved. Third, owing to 
condition ((T^), the limit of Eq. (P7)) as a ^ oo is the electrostatic energy of a 
lattice of point-like charges plus a uniform neutralizing charge within each cell 
(the cubic symmetry being specified by Vq and the a^s). Hence, the limit value 
of the rhs of (P7)) for a ^ oo is equal to Mdi^a /'^a^rs where the M^i^aS are the 
Madelung constants defined below Eq. ([3]). This remark shows that the SLQC 
approximation yields the classical crystal energy in the infinite dilution limit 
because later, after illustrating the results shown in Fig. 2, it will appear clear 
that the limit —^ oo implies that a oo. 

Before concluding the section, we further elaborate the above expressions in 
terms of new integration variable ^ and dimensionless function (f)o{^,a). The 
new quantities are respectively defined by 



and 



r = ^6ai/2. 



M^a) = (8/«o)"'/'0(r,a) = (8/^;o)~'/'0(^ 6ai/2, a). 

1=1 



(29) 



Recalling Eq. ([5]), the tip of vector ^ is confined to the cubic cell Cq of edge 2 
centred at the origin of an orthogonal Cartesian frame and defined as 

Co = {e1-i<ej < 1, i = 1,2,3}. 

00(^1 Oi) is an even function and, due to Eq. pT|) . it obeys to 

Co 



(30) 



Substituting {vo/8Y^^<j)o{^,a) for 0(r, a) in Eq.s (|20p and (|27[) and passing to 
new integration variable ^, one finds that 



(31) 



where the dimensionless quantities T^{a) and Va-{a) are 



Vcr{a) 



]:l=l 



Co 



(32) 



Xa- 



Co 



d,(e-e + 2m) 



Co 



^,0o^ (g,«)/4-l/64 ' 
dA^-i') 



(33) 
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In these two equations, similarly to Mdi,a's definition, index a ranges over 
{1, 2, 3} respectively associated to lattices sc, bcc and fee. Moreover, the partial 
derivatives present in Eq. (|32p refer to variable ^ while the remaining symbols 
are defined as follows: 

diiO^lt^f (34) 

for the sc case, 



K2 



A^,, = 2, and ^(2),, = 1 if , ^l, X2 - (-|^) 



3 



j=l l<J<i<3 

for the bcc case, and 

-3=(^)''', ^(^-^ = 3 and ^(3),, = -1 if ,V ^, X.-{^)"\ 

]=l 1<]<1<3 

for the fee case. 

Numerically it is more convenient to express Uo-(a) in terms of Fourier transforms 
(FT) as we already showed in a first presentation[2T] of the SLCQS approach 
based on quantum field theory but restricted to the only sc case. After denoting 

the FT of (Ci Q;) by 0o^(q, ct), i.e. 

^(q, a) = I d^^e-'^< <j,\ {I a) , (37) 

in appendix A we show that v„{a) can be written as 

Vaia) = wp,cr(a) + VN,cria), (38) 



with 





^/ y-'(0o^(7i-m,a))2 
d2^(m) 




VN.aia) 


= xJ dHf d^e- 

J Co -'Co 


daii-i') 


dhici) 


2 1 2 1 2 
= gi +92 +93 , 






= q • q + 9192 + 9293 + 


9391, 


d^M) 


= q • q- 2(gig2 + 9293 


+ 939i)/3, 



(39) 

(40) 

(41) 
(42) 
(43) 
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and 

Xi' = xi/27r, X2' = \/3x2/47r and xs' = Z'/'xs/eTr. (44) 

After choosing a particular real even function <j)Q{£^,a), which obeys Eq. (|30p 
and vanishes with its first and second partial ^-derivatives on the boundary of 
Co, the best SLQC wave-function is obtained as follows. First, one evaluates 
quantities Ta{a), vp^cr{ct) and VN,cr{a) by Eq.s ([5^ . ([55]) and over a grid 
of a values: ai, . . . ,0;^. Second, chosen an value, by Eq.s (|^. ([5^ . 
and (PO)) one obtains the set of values e^irs, ai) with i = 1, . . . , L. The smallest 
of these values will correspond to a particular i denoted by i. Then, ecr(rs, a^) 
approximates the energy of the fundamental state for the considered SLQC 
wave-function with the considered crystalline symmetry, while represents the 
value of aa{rs). Finally, the value of tr that yields the smallest energy at a fixed 
value determines the crystalline symmetry of the jellium at the considered 
density. 



3 Numerical results 

We report now a numerical illustration of the procedure just described. To this 
aim we define function 4>o{^, a) as follows 



with 



G(C,a) 
C(a) 



EE C(a)e-"«'/(i-«'\ 
= (V^vl/(i/2,0;2a) 



-1/2 



(45) 

(46) 
(47) 



where ^'(1/2, 0; 2a) is a particular value of the confluent Hypergeometric func- 
tion ^'(a, c; z) defined in §6.5 of Ref.[22]. The reported C{a) expression ensures 
that 4>o{^, a) obeys condition ([30|) . Besides, this function as well as all its partial 
derivatives vanish as ^ approaches the border of Co . The factorized expression 
of (I>q(£ a) further simphfies Eq.s ^7^, ([Ml) and (gOl). In fact, as shown in 
appendix B, To- (a) becomes 



Ta{a) ~ k' aT{a) = k', 

with k'i = 3ki, K'2 = 6^2, k! z = 9k3 and 



2 / [d^G{i,a)Yd^ 




T(a) 



M/(-f, -2; 2a) 
4a ^'(i, 0; 2a)' 



io^(q, a) becomes 



6o'(q,a) = nG2(g„a) 



(48) 



(49) 



(50) 
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with 

r-l rl 



G^{q,a)= f d^e-"'^G^{C,a) =2 f d^cos{q^)G^{^,a). (51) 
Besides, after putting for < 77 < 2 

r(-77, a) = r(ry, G^{f] + ^' , a)G^{^', a)dC, (52) 



Eq.s and PO)) respectively become 

,^. n-^i(G^(^m,,c.))^ 

vp,^(a) = X<t2^^^ , , , (53) 

...(«) = -xJ d^rfl^^, (54) 
J2C0 '^'t(?7) 



where 2Co denotes the cubic cell {r/j —2 < rji < 2, i — 1,2,3}, and with the 

d2^(m)s and d^iffjs defined by Eqs. (gllj-jlS]) and (|34| -([36 |l . 
In this way, we must numerically evaluate the ID integral present on the rhs 
of (|48p to determine To- (a), the ID FT defined by Eq. (jSTj) at a set of values 
q = TTm with m = 0, . . . , M to subsequently evaluate series ()53p truncated at 
M, and, finally, r(77, a) at a set of points 771, ... , rjMi lying within the interval 
[0,2] to evaluate 3D integral (|54|) . (By so doing, we use the eveness of r(7y,a) 
and G^{q,a) with respect to rj and q.) The calculations must performed over 
a grid of a values: < ai < . . . < a^. After numerically determining T^{ai) 
by Eq. ([HI), wp,<j(ai) by Eq. (HH]) and UAr^<T(ai) by Eq. ([Ml), we look for the 
minimum value of 

ecr(rs,ai) = 5 1 (55) 

at a fixed value as i ranges from 1 to L. Then, as already reported at the end 
of the previous section, if the minimum occurs at i = i, the SLQC fundamental 
state corresponds to Eq. ([6]) with a = aj and the corresponding energy per par- 
ticle, in Rydberg units, is £a{rs,ct:i). Changing the rg value, i will also change 
and the set of the aj values numerically determines a^irs), i.e. the way 
changes with is determined by the set of values aj relevant to the considered 
lattice symmetry. 

Before illustrating our numerical results, we give some details about numerical 
computations. We considered an a grid that spans the interval [0.01, 500] at 
integer multiples of 0.05 up to a = 0.5, of 0.1 up to a = 2, of 1 up to a = 20, 
of 10 up to a = 100 and of 100 up to a = 500. The numerical accuracy of the 
results, depending on the values of AI, Mj and the number of points used to 
evaluate quantities T{a), G^{q,a), T(ri,a) and Ucr(a), was tested at a = 0.01 
and 100. The first six significant digits of T(a) do not change passing from 10^ 
to 10^ points. The first 5 digits of the FTs evaluated with 10^ or 10^ points do 
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not differ and agree with their asymptotic formula given by Eq. (96) of [21] [note 
that here the correct argument of the sine function is (2 + \/2qa — tt/S)]. The 
first five digits of the resulting up.(j(a) values, obtained truncating the series at 
M = 50 or M = 100, also coincide. The calculation of VM,cT{<y) requires more 
care. The evaluation of r(77) with integration steps of 10~^ or 10~^ leaves the 
first five digits unchanged. However, the evaluation of the remaining 3D integral 
over fj cannot be done with an integration step A?7 smaller than 10~^ because 
the required CPU time becomes very large. Thus we made three runs with Ary 
respectively equal to 1/500, 1/1000 and 1/1500. The first four digits do not 
change even at large as, i.e. a > 100, which is the critical region. The ratio 
of the resulting errors, defined as (second -fist)run divided by (third - first)run, 
was found fairly equal to 3/2. In this way we numerically extrapolated the re- 
sults to A77 0. The resulting VN_a{a) values ought to have five correct digits. 
vp,aici) and VN,a-ic() always have opposite signs. At large as, they have the first 
two digits equal. Hence, the final Va{a) values have the first three digits exact 
for all the considered values of a. 

We pass now to illustrate our results. These are shown in Fig.s la, lb and 2. 
Fig.s la and lb respectively show the energy per particle (in Rydberg units) 
versus Ts, in the high density region, and vs Lgio{rs) in the low density one. 
In fact, the continuous and dotted curves refer to Bloch's expression ^ and 
Carr's expression ([5|), respectively. The long-dash curves refer to the optimized 
SLQC wave-function with symmetry bcc, and the short-dash curves to the sim- 
ple cubic optimized SLQC wave-function. The full triangles and circles are the 
values respectively obtained for the fluid and bcc crystalline phases by QMC 
calculations [TTl \T2\ [T5l IT6] . (We have not reported the curve relevant to the 
fee SLQCS solution in order to not overcrowd the figure. It lies close to the 
bcc solution.) The short-dash curve crosses the continuous one at Vg = 28 and 
for greater values it lies below the continuous curve. On the basis of the 
Ritz-Rayleigh principle one concludes that the SLQCS with the sc symmetry is 
closer to the true fundamental state than \Fp). In other words, as the density 
decreases, the jellium passes from the polarized fluid phase to the sc crystalline 
one at Ts = 28. The figures also show that long-dash curve passes from above to 
below the continuum one at = 38. This means that in the region 28 < < 38, 
the sc phase is more stable than the fluid which in turn is more stable than the 
bcc. In the region > 38, the bcc phase is more stable than the fluid and is 
less stable than the sc up to = 500. The figures make also evident a property 
that appears to have been overlooked by most textbooks: Carr's approximation 
appears to be surprisingly accurate throughout the full density range since it 
fairly agrees with QMC results even at high densities. 

We pass now to the illustration of figure 2. The full triangles and the close dotted 
line respectively represent a2{rs) and ai{rs) on a log- log scale. The apparent 
step behaviour is an artifact of the chosen a-grid whose values are not equally 
spaced. One should note, at large LgiQ(rs)s, the approximate linear behaviour 
of L(7io(a) with a slope corresponding to have oc Vs^^^. The continuous and 
the long-dash curve reported on the top of the figure respectively plot U2(a) and 
fi(a) vs. Lgio(a). These curves can also be considered as the plots of (a (rs)) 
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Figure 2: The dotted step-line and the fuU triangles are the plots of 10 x 
Lgw{aa) vs. Lgio{rs) for the so and bee case; the continuous and the long-dash 
curves (convex and monotonically decreasing) those of r^^^[eairs)+Mdi,a/'r's] vs. 
Lgw{rs) for the bcc (cr = 2) and sc case (cr = 1) [the corresponding vertical scale 
is on the left and the energy units are ryd]; the continuum and the long-dash 
sigmoidal curves those of V2{ce) and vi{a) vs. Lgio{a) (the relevant horizontal 
and vertical scales are the top and the right ones); finally, the parabolic dotted 
curve that of Lgio{T2{a)) vs. Lgio{a) (its vertical scale is the left one). 
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vs. Ts for the noted property that (x r^^^'^. In this way, the figure makes 
it evident that the SLQCS Va{rs)s approach the relevant Madelung values as 
Ts — > oo. The remaining two monotonically decreasing convex curves are the 

plots of r^^'^[£a-{rs) + Mdi.a/fs] em vs. for the bcc (continuous curve) and 
the sc (dotted) symmetry. The constant behaviour observed at large values 
indicates that 

eArs)~ — + -17I+o(-^), (56) 

the Co-s being appropriate constants that from the figure appear to be equal 
to 3.6 and 2.9 for cr = 2 and 1, respectively. In appendix B we show that 
the SLQCS procedure always yields an energy per particle that asymptotically 
behaves as reported in Eq. at large Vg. It is stressed that this behaviour 
coincides with that of Eq.s ([U and (O, even though the numerical coefficients of 

3/2 

the Ts contribution are different. In this respect, we recall that the crystalline 
HF solutions investigated in Ref.[10j also behave as in Eq. ([56| and from Fig. 4 
of this paper it appears that C2 ~ 3. We once more underline that contribution 
—Mdi^2li's comes naturally out by the HF and the SLQCS procedure [that is 
equivalent to a Hartree or a HF equation as explained below Eq. while in 

Wigner's and Carr's derivation of Eq.s ([4]) and ([5]) it was put there because it 
is the energy of the 'static' unperturbed Hamiltonian. 

Finally, the dotted parabolic curve plots L(7io(t2(q!)) vs LgiQ^a). It shows a 
linear behaviour in the outermost a range that follows from Eq. (|49p as we 
explicitly show in appendix B. 



4 Conclusion 

As Wigner first pointed out, an overall neutral one component plasma of elec- 
trons shows, at T = O^K, the interesting feature of being in a crystalline phase as 
the particle number density becomes smaller than a particular value. This fea- 
ture can nowadays be considered well assessed because it was confirmed by HF 
and QMC calculations (though the transition density is not accurately known 
yet). The fact that quantum Coulombic crystals exist at very high dilution, 
on physical grounds, has the important consequence that their behaviours must 
approach those of the corresponding classic Coulombic crystals. In fact, in the 
infinite dilution limit, the inter-electron mean distance becomes infinitely large 
with the consequence that: the overlapping among the the wave- functions of 
different electrons is expected to vanish and electrons become distinguishable 
in the sense that one can speak of an electron that occupies a well definite 
region/cell of the crystal. The SLQCS approach exactly captures the last two 
features. In fact, the SLQCS is completely antisymmetric but the support prop- 
erties of the involved wave- functions make the fully antisymmetric wave- function 
equivalent to the Hartree one. Consequently, the electrons behave as particles 
each of them being confined to a single cell of the crystal. The determination of 
the best SLQCS, equivalent to solve the corresponding Hartree or Hartree-Fock 
equation, proceeds by a simple variational procedure. In fact, we worked out 



16 



the expression, in an integro-difFerential form, of the energy per particle, given 
by Eq. ((3T|) together with Eq.s ((32|) - ((43|) . in terms of the basic function (j){r, a) 
different from zero only within a single primitive cell and depending on a real 
parameter a. The solution was obtained by looking for the minimum of the en- 
ergy with respect to a at each fixed value. The results illustrated in the three 
figures show that, at large r^s, the leading term of the SLQCS energy per parti- 
cle coincides with that of Eq.sQ and (O with an important difference: Wigner 
and Carr took this contribution from the classical Coulombic crystal's value 
while the SLQCS procedure directly produces this term. Moreover, as shown in 
appendix B, this procedure also implies that the next to the leading asymptotic 
term of e<T(?'s) behaves as Cr,lr^l'^ . The dependence coincides with that of 
Wigner's and Carr's formulae but the numerical values of coefficient Ca are dif- 
ferent. These values are also different from the values obtained solving the HE 
equation|10j. Hence, in the very low density region, the differences among the 
Wigner/Carr, the SLQCS and the HE approximations of the jellium model set 
only into at the level of the 0(rs~^^'^) term. We consider this result as the main 
conclusion of this analysis. The implications of this conclusion are: Wigner or 
Carr formulae are equally accurate in the far Vg region though the dominat- 
ing contribution was not quantum-mechanically derived; the SLQCS approach 
quantum-mechanically derives this contribution and therefore it is there as ac- 
curate as the HE equation but of simpler application; as the system becomes 
denser the numerical differences in coefficient Ca make the approximations no 
longer equivalent. One expects that overlapping effects be no longer negligible 
and that the SLQCS approach drastically deteriorates in comparison to the HE 
one. 



A Conversion of Eq. (jSl) in Eq.s ([MD and ([40]) 



To write Eq.s (|33p in a form more convenient for numerical computations, we 
first introduce function So(C)j equal to 1 if ^ lies within Co and equal to zero 
elsewhere. After putting 

^{l a) = 0o'(^, a) - eo(0/8, (57) 

one easily proves that the quantity inside the square brackets on the rhs of 
Eq. ([55)1 can be written as 

y / dH f rf3^-,M;^M>i _ [ f ^z^.'^igs^liLA 

i^Jco J Co d^(|'-^' + 2m) J Co J Co d„{£^-S!') 

(58) 

where the sum includes now the term with m = 0. Since (jioi^jO:) and 6o(C) 
are equal to zero outside Co, integrals in Eq. ([55]) can be extended to the full 
R^. The ET transform of \E'(^, a) wiU be denoted by ^(q, a). Both 0o(C, a) and 
0o(C) are real and even function of ^. Then \I'(q, a) is an even function of q. 
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Using the FTs, the generic term of the series can be written as 



V'T^) Jr3 Jm - f + 2m) 7 J 

The change of the integration variable ^^77 = ^ — + 2m allows us to perform 
the ^' integration so as to convert the previous expression into 



a 



The integral over fj yields 



1 /" ,3^ 1 47r uj„ 



(2^)3 'd„{rf) (2^)3rf2^(q) 

with 



= 1, c^2 = V3/2, c^3 = 2^/3/3, (60) 



and the (Pa{f^)^ defined by Eq.s ((1T1) - (|^51) . In this way, the series present in 
Eq. ([55)1 converts in 

(2^)3 J dKici) 
Using the mathematical identity [S] 

E^^''"" = ^E^(q-^™) (62) 

m 111 

the integrals over q in Eq. (I6T|) are immediately evaluated if m ^ 0, while the 
contribution relevant to m = is equal to zero because the resulting integrand 
\I''^(q, a)/rf^cr(q) vanishes at q = 0. In fact, condition ((30|) and the definition of 
&o{0 imply that 

^{0, a)= f d^^ [0o'(e, «) - 60(6/8] = 0. 

J Co 

Moreover the eveness and the reality of 00 (Cj «) and 60 (0 ensure that 

^(q, a) = 0o'(q, «) - e'o(q)/8 « 0(q • q) 
at small |q|. This property implies that 
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and one concludes that no contribution to the sum over m arises from the term 
with m = 0. Thus, the series present in Eq. ((58)) is equal to 

u^yn-m,a)^ (63) 

A further simplification follows from the analytic expression of 0o(q). This 
reads 

3 ^1 3 . / ^ 

sm(gj) 



0o(q) = n / e-'*«d^ = 8n 



.=1 

Since 0o(7rm) = if m 7^ 0, we can replace ^^(7rm, a) with (0o^(7rni, a))^ in 
Eq. dMl) and finally write Eq. ([58]) as 



'3c f j3^5'/'o^(^,a)</>o^(C',a) 



d'q dY -" vX j • (64) 

Co J Co a<T(?-C) 



Expression p8|) for UCT(a) and Eq.s ([39]) and (|40|) for up,cr(a) and UAr,o-(a) imme- 
diately follow from Eq. (|64p recalling that Eq. ([58]) is the content of the square 
brackets in Eq. 



B Asymptotic behaviour of Cairs) at large 

First of all Eq.s (|47| and ([49]) are derived as follows. The condition that (/'o(Ci ck) 
be normalized requires that G(^, a) be normalized and this implies that 

C-^a) = 2 r e-2"-'/(i-"')dx. (65) 



With the change of the integration variable: x y^/ -^1 + y, the above expres- 
sion converts to 

/•oo 

C-\a)^ 2;-i/2(l + y)-3/2e-2"!'dy. (66) 







Recalling the general definition of the confluent hypergeometric function '5 (a. 6; z) 

1 

^{a,b;z)=—— e-^'H''-^ {1 + ty-''-Ut, (67) 



r(a) Jo 

reported in §6.5 of Ref. [22\ . from Eq. (|66p immediately follows that 

C-2(a) = V^^i^, 0; 2a), (68) 
which is equivalent to Eq. (|47p . In the same way, 

2/ [a.G(x,«)]'d^ = 2C2(a) / [^^^^^ ^ fdx (69) 
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and, by the previous change of the integration variable, one finds 



2C\a) / [2a2e-2-yyi/2(n_ 2^)3/2] 



= C2(a)V^vl/(_ 2,2a) 



-2, -2a) 



4a*(i, 0, 2a) 



(70) 



i.e. Eq. (j49| . 

Finally we show that the SLQCS approximation implies that tail's) at large r, 
behaves according to Eq. (|56p . Once we have determined acr{rs), from Eq. ([3T 
follows that 



2ao 



r(a(rs)) t7CT(a(?'s)) 



while the aCT(''s)s are determined solving the equations 

d ( Tcr(a) ■Llcr(a) ' 



0. 



(71) 



(72) 



Illustrating Fig. 2's results, we already noted that limQ^oo ^^{0) — —Mdi^a and 

1 /2 

that a Ts at large Vg. Besides, based on the fact that the Va{rs)s approach 
their limit values from the above, after putting (3 = 1/a it appears reasonable 
to assume for the Ucr(a)s, as /3 0+, the following asymptotic behaviour 



(73) 



(74) 



Ccr being a positive constant. The asymptotic behaviour of To- (a) at large a is 
easily obtained from that of ^(a, ,6; 2a), reported in §6.13.1 of Ref.[55]. One 
finds that 

3 9/3 
_/3 ^ 2 ^ 16 

Eq.s ([75)1 and ([71)1 allow us to evaluate the derivatives present in Eq. ([7^ so 
that this equation, to the leading order, converts into — Ko./(/3^r^) + Ca/vs — 0. 
The solution is 

a = \J CaTs/Kcr. (75) 

Its substitution in Eq. ([71)) yields the first two leading terms of the asymptotic 
expansion of ea{rs) at large r^, ie 



e<T(rs) 



e 

2ao 



3/2 



(76) 



This result shows that the SLQCS approach implies that, at very high dilu- 
tion, the leading term of the energy per particle is the classical Madelung value 

3/2 

and that the leading correction to this term decreases as with a positive 

numerical factor equal to ^/K^c^. Fig. 2 results indicate that this factor nearly 
equals 2.9 and 3.6 in the sc and bcc case, respectively. 
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